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First-principles calculation of the thermal properties of silver 
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The thermal properties of silver are calculated within the quasi-harmonic approximation, by using 
phonon dispersions from density-functional perturbation theory, and the pseudopotential plane- wave 
method. The resulting free energy provides predictions for the temperature dependence of various 
quantities such as the equilibrium lattice parameter, the bulk modulus, and the heat capacity. Our 
results for the thermal properties are in good agreement with available experimental data in a wide 
range of temperatures. As a by-product, we calculate phonon frequency and Griineisen parameter 
dispersion curves which are also in good agreement with experiment. 
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I. INTRODUCTION 

The study of the temperature dependence of the prop- 
erties of materials requires a proper account of nuclear 
motions. Within the framework of density functional 
theory (DFT) 0, a major breakthrough in this field has 
been opened by the introduction of ab-initio molecular 
dynamics by Car and Parrinello Far from the melt- 
ing point, however, a more conventional (and yet largely 
unexplored in practice) approach based on lattice dynam- 
ics proves to be both more accurate and computationally 
efficient. In the harmonic approximation, the crystal free 
energy is calculated by adding a static contribution — 
which is accessible to standard DFT calculations — to a 
dynamical contribution which is approximated by the free 
energy of system of harmonic oscillators corresponding 
to the crystal vibrational modes (phonons). The latter 
is nowadays conveniently calculated by using density- 
functional perturbation theory (DFPT) In the 



JtJ, some anharmonic 
Towing phonon frequen- 



quasi-harmonic approximation 
effects can be accounted for by af 
cies to depend on crystal volume. Among other advan- 
tages, the (quasi-) harmonic approximation allows an 
explicit account of quantum effects on nuclear motion, 
which can be important below the Debye temperature. 
Furthermore, analysis of the normal vibration modes 
and of their individual contribution to the free energy 
can explicitly reveal the mechanism driving the thermal 
expansion, phase transitions, and the crystal stability. 



The main concern about the lattice dynamics method is 
the range of validity of the quasi-harmonic approxima- 
tion. Calculations based on various semi-empirical mod- 
els P-pd| as well as on first-principles methods [jl 12 - T^] 
demonstrate that the quasi-harmonic approximation pro- 
vides a reasonable description of the dynamic properties 
of many bulk materials below the melting point. Very 
recently, first-principles calculations on the thermal ex- 
pansion of some simple s-p metals indicate that the treat- 
ment of anharmonic effects at the quasi-harmonic level 
provides a remarkable good description of the structural 
and elastic properties of these materials up to their melt- 
ing points Q . 

In this paper we apply the quasi-harmonic approxima- 
tion to the study of the thermal properties of the 4d noble 
metal Ag, such as thermal expansion, heat capacity, and 
temperature dependence of the bulk modulus. To this 
end, we first calculate the phonon dispersion curves as 
functions of volume, by using DFPT ||. Our results 
demonstrate that all these quantities can be accurately 
predicted from the present parameter-free method in a 
wide range of temperatures. 

The paper is organized as follows. In Sec. II, we briefly 
outline our computational framework, as well as some 
definitions concerning the physical quantities we have in- 
vestigated. The results of our calculations are then pre- 
sented and discussed in Sec. III. Finally, in Sec. IV we 
give our conclusions. 
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II. THEORY 

A. Equation of state and thermal expansion 

For a given temperature, T, and volume, V, the equi- 
librium state of an extended system such as a crystal is 
determined by the condition that the Helmoltz free en- 
ergy, 



F(V, T) = U - TS, 



(1) 



is at a minimum with respect to variations of all possible 
internal degrees of freedom, such as e.g. atomic arrange- 
ments and electronic states. Here U and S indicate the 
internal energy and entropy, respectively. The equation 
of state of the system is obtained from Eq.( |l|) by equat- 
ing the pressure to minus the volume derivative of the 
free energy: 



P = 



'dF\ 
K dV) T 



(2) 



In the quasi-harmonic approximation, F is given by 
F(V,T) = E(V) + F vib (io,T) 



(3) 



where E is the static contribution to the internal 
energy — which is easily accessible to standard DFT cal- 
culations, F V ib represents the vibrational contribution to 
the free energy, and u>j(qi) is the frequency of the j-th 
phonon mode at wave vector q in the Brillouin zone (BZ). 
Anharmonicity is explicitly, though approximately, ac- 
counted for by allowing E(V) to deviate from a quadratic 
behavior and by letting the phonon frequencies depend 
on volume. Since the temperatures considered here are 
well below the electronic energy scale, the contribution 
of the electronic excitations to the thermal expansion is 
negligible and is not included in the present work. The 
equation of state (|3|) can now be written in the form 



dE 



dF v 



ib 



7.7 (q) 



dV iv j (q) ' 



(5) 



and £(ojj(q)) is the mean vibrational energy of the (q, j) 
phonon given by 



£( w .j(q)) = ftwj(q) 



i 



i 



2 exp(7k^(q)/fc B T) - 1 



(6) 



The thermal expansion is obtained directly from the 
equation of state (0) and the volume thermal expansion 
coefficient is defined as 



ay 



V [df 



(7) 



The temperature dependence of the bulk modulus is ob- 
tained from 



B(T) 



v 



V dV 2 ) , 



d 2 E 
dV 2 



V 



d 2 F mb (u J ,T) 
dV 2 



(8) 



Due to anharmonicity, the heat capacity at constant 
pressure, C p , is different from the heat capacity at con- 
stant volume, Cy. The former, which is what experi- 
ments determine directly, is proportional to T at high 
temperature, while the latter goes to a constant which 
is given by the classical equipartition law: Cy ~ 3Nk%, 
where N is the number of atoms in the system. The 



relation between C p and Cy is 



C„ — Cy — 



r(T)BVT, 



and Cy is given by 



q j 



V 2k B T J smh 2 (huj j {q)/2kBT) 



(9) 



(10) 



B. Computational details 



dV 



| + I^^ 7j(q )%(q)), 



where jj (q) is the Gruneisen parameter corresponding to 
the (q, j) phonon mode, defined as: 



The static total energy, E(V), and phonon frequencies, 
(4) w j(q) 5 are calculated by using DFT and DFPT respec- 
tively, within the local density approximation fFfll. We 
use separable norm-conserving pseudopotentials JL8 19 
together with a plane-wave basis set up to a kinetic- 
energy cutoff of 55 Ry. Sums over occupied electronic 



2 



states are performed by the Gaussian-smearing special- 
point technique 20 23], using 60 k points in the ir- 
reducible wedge of BZ. Phonon frequencies are calcu- 
lated on a (444) regular mesh and Fourier-interpolated 
in-between. This Fourier interpolation amounts to in- 
cluding real-space inter-atomic force constants up to the 
nincth shell of neighbors. 



III. RESULTS 

Figure [j] shows the static total energy per atom, E(V), 
as a function of the lattice constant a — (AV) 3 (Ag is a 
face-centered cubic metal). Our data are fitted to a Mur- 
naghan's equation of state [ p2[ . The resulting lattice con- 
stant, ao — 4.05 A, bulk modulus, B = 1.28 Mbar, and 
pressure derivative of the bulk modulus, dB/dp = 5.66, 
agree well with previous theoretical calculations |^3| . For 
comparison, the room-temperature experimental data 
are: a = 4.08 A fH and B = 1.01 Mbar @. 

According to Eq. (|4|), in order to obtain the equation 
of state one must first calculate the phonon band struc- 
ture as a function of volume. In Figure ^| we display 
the phonon dispersion curves as calculated along several 
symmetry directions at the minimum of the static en- 
ergy. Experimental data at room temperature [^5| are 
reported for reference. The effect of temperature on the 
phonon dispersions will be discussed later in this paper. 
Figure || shows the calculated dispersion curves of the 
mode Griineisen parameters of Ag, as defined by Eq. (||), 
along the same symmetry directions. The dispersions are 
discontinuous at the BZ center as a consequence of the 
anisotropy and polarization dependence of the sound ve- 
locities. The Griineisen parameters of silver are positive 
throughout the BZ for all branches, thus implying that 
there is no anomalous negative thermal expansion at low 
temperature, as in Si p6| . The averaged Griineisen pa- 
rameter of silver is 2.6, in agreement with the experiment 
value of 2.5 §7|. 

With the static total energy and the Griineisen param- 
eters in our hands, we can set up the equation of state via 
Eq. (||). Figure H shows the pressure p(V, T) as a function 
of the lattice parameter a for several temperatures. It can 
be seen that for a given pressure, the lattice constant in- 
creases with the temperature. At room temperature and 
zero external pressure, the calculated lattice constant is 
ao = 4.07 A, which is closer to the room temperature 
experimental value of 4.08 A [^4| than the result of 4.05 
A derived from the static total energy. As the tempera- 



ture is increased, a critical temperature T m exists above 
which p(V, T) no longer intersects the p = line and, 
within the quasiharmonic approximation, the crystal be- 
comes mechanicaly unstable due to the vanishing of the 
isothermal bulk modulus — V(dp/dV)T at T m . In Silver 
this lattice instability occurs at 1370 K, not far from the 
experimental melting temperature (1234 K) jlSL 

The thermal expansion can be derived directly from 
the equation of state. We have calculated the linear ther- 
mal expansion which is defined El] as 



Aa _ a (T) - a {T = 293) 



«o 



a (T = 293) 



(11) 



The results are shown in Figure ||. The agreement be- 
tween the theoretical and experimental results is very 
good. Although the calculated equilibrium lattice con- 
stant is slightly different from experimental measure- 
ment, the temperature dependence of the relative vol- 
ume changes are described accurately by the present 
method. The good agreement between theory and ex- 
periment holds not only at low temperature, but also 
near the melting point. 

Given the volume as a function of temperature, the 
temperature dependence of phonon frequencies can be 
approximately estimated through their volume depen- 

' j(T)-w(T=293)) 



dence: 



w(T=293) 



37s, where 7 is the Griineisen 



parameter ( Eq .( |5|)) and e is the linear thermal expan- 
sion, (Eq.( |ll|)). In the present case, given the fact 
that e(T = 0) ~ —5 x 10~ 3 , and that Griineisen pa- 
rameters are tipically 7 w 2 ~ 3, phonon frequencies 
at room temperature are lower than their T — value 
by 3 ~ 5%. Figure || shows the theoretical phonon dis- 
persions at room temperature (full lines) together with 
the dispersions, which are already shown in Figure |^ and 
calculated at the static equilibrium (dashed lines in Fig- 
ure ^) . The comparison with experimental data at room 
temperature (open circles) p5| is clearly improved and 
now rather satisfactory. 

Figure |t] shows the temperature dependence of the bulk 
modulus as calculated from Eq.( ^). Bo denotes the bulk 
modulus obtained from the static total energy by ne- 
glecting the lattice vibrations. At room temperature, we 
obtain a bulk modulus B = 1.16 Mbar. Comparing this 
value with the result derived from static LDA calcula- 
tions (Bo = 1.28 Mbar), one sees that the agreement 
with the experimental result (B = 1.01 Mbar) [Jl6| is 
significantly improved. 

The calculated heat capacity C p and Cy are shown in 
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Figure |[ It can be seen that below the Debye tempera- 
ture (9d = 215 K for silver |p7j ) the difference between 
C p and Cy is very small, while at high temperature, the 
heat capacity at constant volume Cy approaches to the 
classical value 2.49 Jmol _1 K , while the heat capacity 
at the constant pressure increases monotonously with the 
temperature. The available experimental data [ p9[ for C p 
are shown as circles. The agreement between theory and 
experiment is remarkable in a wide range of temperatures 
also in this case. 



IV. CONCLUSIONS 

In the present paper, we have calculated the thermal 
properties of silver, such as thermal expansion coefficient, 
Griineisen parameters, bulk modulus, and heat capac- 
ity, using the the quasi-harmonic approximation within 
density-functional theory. The equilibrium lattice con- 
stant is obtained from the equation of state constructed 
by the free energy. The volume dependence of frequen- 
cies is calculated from the density-functional perturba- 
tion theory. The obtained results for the investigated 
thermodynamic quantities are in good agreement with 
the available experiment measurements. The calculation 
suggests that the anharmonic properties of silver can be 
accurately calculated from this first-principles approach 
in a wide range of temperature. The application of the 
present method to the study of other properties of ma- 
terials such as high pressure effect, crystal stability and 
phase transitions, is straightforward. 
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FIG. 1. Static total energy per atom as a function of lattice 
parameter a 
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FIG. 2. Calculated phonon dispersion curves at the lattice 
parameter corresponding to static equilibrium. Experimental 
neutron-scattering data |^5[ are denoted by circles. T and 
L represent transverse modes and longitudinal modes respec- 
tively. 



FIG. 3. Calculated dispersion curves of the mode 
Griieisnen parameter 7j(q) of silver along some symmetry 
lines in BZ. T and L denotes the transverse modes and longi- 
tudinal modes respectively. 
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FIG. 4. Applied pressure P(V, T) for silver as a function 
of lattice parameter a (here V = a 3 /4) for different tempera- 
tures. 



FIG. 5. Temperature dependence of the linear thermal ex- 
pansion for Ag. Solid curve is the calculated result and the 
circles represent experimental data from Ref. |2q]. 
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FIG. 6. Solid lines: calculated phonon dispersion curves FIG. 7. Calculated temperature dependence of the ratio of 

of silver at T — 293 K. Dashed lines: calculated phonon bulk modulus B(T) to Bq (obtained from static total energy) 
dispersions at the lattice parameter corresponding to static for silver, 
equilibrium. Experimental neutron-scattering data |^5j are 
denoted by circles. T and L represent transverse modes and 
longitudinal modes respectively. 
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FIG. 8. Calculated temperature dependence of heat ca- 
pacity of Ag at constant pressure [C p ) and at constant vol- 
ume (CV). The experimental data for C p [^] are denoted by 
circles. 
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